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1  Introduction 


Radar  systems  engineers  often  obtain  a  taper  for  a  planar  array  as  the  separable  product  of  one¬ 
dimensional  tapers  in  two  orthogonal  directions.  While  this  is  rarely  appropriate  for  final  design 
work,  it  is  often  adequate  for  various  preliminary  designs,  for  example  for  simulations  or  other 
calculations  where  obtaining  rough  array-performance  numbers  will  settle  some  matter  at  hand.  Of 
course  the  analogous  separable-filter  idea  is  well  known  in  image  processing. 

In  the  image-processing  community,  however,  a  second  design  procedure,  the  McClellan  trans¬ 
formation,  for  creating  a  two-dimensional  filter  from  a  one-dimensional  prototype  filter  is  both  well 
known  and,  by  some  authors’  claims,  widely  used.  The  transformation  from  one  dimension  to  two 
is  mediated  by  a  transformation  kernel — a  “spreading  function”  in  this  report — whose  characteris¬ 
tics  affect  the  final  two-dimensional  design  as  much  as  does  the  one-dimensional  prototype.  The 
extra  degrees  of  design  freedom  represented  by  the  choice  of  transformation  kernel  represent  flex¬ 
ibility  not  available  in  the  simpler  separable-filter  approach,  so  it  is  natural  to  wonder  whether  the 
McClellan  transformation  might  be  useful  in  the  design  of  planar  array  tapers,  either  as  an  alterna¬ 
tive  to  separable  tapers  for  quick-and-dirty  designs  or  perhaps,  if  adequate  performance  turns  out 
to  be  available,  in  final  designs  as  well.  Remarkably,  a  search  of  the  literature  for  work  on  using 
the  McClellan  transformation  to  design  array  tapers  turns  up  nothing  but  a  short  letter  by  McNa¬ 
mara  and  Botha  HI  pointing  out  that  the  approach  can  be  used  to  design  wide-coverage  beams  for 
satellite-communication  arrays. 

The  purpose  of  this  report  then  is  to  present  the  McClellan  transformation  somewhat  tutorially, 
to  examine  it  in  an  array  context,  and  to  explore  the  associated  array-taper  design  space  on  a  pre¬ 
liminary  basis.  It  must  be  cautioned  though  that  no  attempt  is  made  here  to  systematically  explore 
that  design  space,  to  do  careful  trade  studies,  or  to  come  up  with  precise  numbers  to  quantify  the 
array-design  capabilities  of  the  approach.  Though  appropriate  performance  numbers  are  provided 
for  all  design  examples  presented,  the  intent  here  is  really  qualitative  exploration  of  enough  parts 
of  the  space  to  perhaps  suggest  to  future  researchers  what  systematic  explorations  might  ultimately 
be  appropriate,  and  this  approach  yields  design  experiments  that  stop  well  short  of  final  answers. 
This  is  intentional,  as  detailed  quantitative  trade  studies  are  better  performed  in  specific  application 
contexts. 

No  formal  distinction  is  made  in  this  report  between  transmit  and  receive  applications,  but  most 
of  the  design  examples  presented  are  far  more  natural  in  a  receive  context.  Those  few  examples 
explicitly  addressing  beam  broadening  would  of  course  be  of  more  interest  in  a  transmit  context, 
though  only  the  still  unusual  type  in  which  amplitude  tapers  are  permitted.  Phase-only  transmit 
tapers  are  of  no  interest  here,  because  a  phase-only  weight  characteristic  is  not  preserved  by  the 
McClellan  transformation:  a  phase-only  prototype  transforms  to  a  two-dimensional  taper  that  lacks 
the  phase-only  property. 

A  straightforward  IEEE  Xplore  subject  search  turns  up  a  rich  engineering  research  literature 
on  the  McClellan  transformation,  only  a  few  papers  from  which  will  be  discussed  explicitly  here. 
The  variant  of  the  transformation  explored  here  was  first  presented  by  McClellan  E)  in  a  paper  that 
remains  a  solid  starting  point,  even  though  its  detailed  design  examples  are  limited  to  the  square 
sampling  lattice  in  two  dimensions.  A  hexagonal- sampling  environment,  such  as  used  in  most 
serious  planar-array  work,  is  dealt  with  briefly  in  the  presentation  of  the  McClellan  transformation  in 
Dudgeon  and  Mersereau’s  classic  text  (3)  on  multidimensional  signal  processing,  but  the  notational 
approach  to  the  hexagonal  lattice  used  is  clumsy  and  not  easy  to  follow.  The  matrix-based  approach 
of  this  report  is  easier  to  work  with. 

One  very  early  but  now-classic  treatment  of  the  McClellan  transformation  is  the  two-part  paper 
by  Mersereau  and  Mecklenbrauker  (and  Quatieri,  in  the  case  of  the  first  part)  00.  Part  I  focused  on 
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the  design  of  the  spreading-function  transformation  kernels,  but  the  year  was  1976,  when  minimiz¬ 
ing  computational  complexity  was  important  in  a  way  that  it  no  longer  is,  and  the  paper’s  resulting 
consideration  of  only  the  very  simplest  examples  makes  it  of  minimal  actual  relevance  today.  Part  II 
worked  through  the  details  of  computational  efficiency  in  realizing  the  transformed  filters,  but  the 
filter  structures  considered  possessed  unfortunate  coefficient-sensitivity  properties  and  so  became 
largely  irrelevant  in  the  following  year  when  McClellan  followed  up  on  his  original  work  by  pub¬ 
lishing  the  far  superior  Chebyshev  structure  on  which  this  report  is  based  0.  Part  II  also  concluded 
a  bit  too  much  from  the  limited  set  of  design  examples  that  had  been  explored  by  that  early  date,  for 
example  asserting  that  “the  technique  is  only  truly  valuable  for  designing  linear  phase  filters  with 
quadrilateral  symmetry.”  This  report’s  final  design  example,  in  Section [7.3. 5[  certainly  shows  that 
to  be  incorrect. 

The  McClellan  transformation  is  traditionally  considered  an  approach  to  2D  filter  design,  but  in 
fact  it  is  just  as  suitable  for  3D  design.  Antialiasing  filtering  for  staggered  sampling  of  planar  arrays 
(6),  for  example,  requires  cone  filters  in  3D,  and  the  McClellan  transformation  has  recently  been 
examined  for  3D  cone  filters  (7). 

In  any  case,  the  reader  is  encouraged  to  consider  the  present  report  to  be  suitable  background 
reading  prior  to  tackling  the  literature,  rather  than  the  other  way  around,  if  only  to  start  out  with  a 
single  notational  approach  rather  than  many. 

The  next  section  is  an  introduction  to  the  key  concepts  in  the  context  of  a  simpler  lD-to-lD 
transformation.  Then  Section  [3j  the  bulk  of  the  report,  addresses  the  design  of  two-dimensional 
array  tapers.  The  Section]?] conclusion  briefly  synopsizes  the  main  lessons  learned  from  this  effort. 
Finally,  readers  wondering  whether  the  McClellan  transformation  can  be  experimented  with  easily 
should  have  a  quick  look  at  the  Fig.[6]matlab  code  on  page[l4]  The  answer  is  yes,  but  only  once  the 
ideas  and  notation  are  understood. 


2  The  idea  and  its  application  to  ID  filters 

2.1  The  Chebyshev  forms  of  a  ID  linear-phase  FIR  filter 

A  ID  zero-phase  filter,  which  is  just  a  ID  linear-phase  filter  with  impulse  response  centered  at  the 
origin,  has  a  transfer  function  of  the  form 


M 

H{z)  =  hQ  +  YJhk{z~k  +  zk) 

k= 1 

or,  equivalently,  an  impulse  response  that  can  be  written 

M 

h(n)  —  hoS(n)  +  hk(8(n+k)  +  8{n— fc)). 

k= 1 

Let  us  rewrite  this  with  the  dependence  on  time  n  suppressed  for  brevity.  We  can  use  Sk  rather  than 
the  usual  S(n  —  k)  to  denote  a  unit-sample  function  at  time  k. 

h  =  h06°  +  Y^  2% - •  (1) 

k= 1 

Here  5 0  and  each  8~k  +  8k )  are  building  blocks  from  which  the  filter  impulse  response  is  created 
using  simple  weight-and-sum  linear  combining.  We  will  be  using  a  variety  of  such  building  blocks 
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Fig.  1  Chebyshev  direct  form  of  a 
linear-phase  FIR  filter,  in  zero-phase 
form.  Causality  shimming  will  be 
needed  at  each 


heavy  vertical  line 


Fig.  2  Reversed  Chebyshev  form  of 
a  linear-phase  FIR  filter,  in  zero-phase 
form.  Causality  shimming  will  be 
needed  at  each 


heavy  vertical  line 


in  the  discussion  to  follow,  and  it  will  prove  convenient  to  have  them  standardized  to  have  unity  DC 
gain.  This  is  the  reason  for  introducing  the  factors  of  two. 

We  can  synthesize  most  of  these  DC-normalized  impulse-pair  building  blocks  needed  from 
closer- spaced  ones.  We  first  convolve  with  +  51  to  spread  a  building  block  out  in  time,  and  we 
then  subtract  unwanted  terms. 

5-2+52  =  (r1+^1)*  ^±^-<5° 


2 

S~k+8k 


=  {8  x+  (51)  * 


£-(&-!)  +  Sk~l  +  Sk~ 2 


(2) 


the  latter  for  any  k  >  2  or  in  fact  for  any  k  >  1  because  the  first  equation  here  is  just  the  k  —  2 
special  case  of  the  second. 

We  can  clean  this  up  with  some  definitions.  We  give  the  spreading  function  +  S1  and  the 
DC-normalized  building  blocks  \  ( S~k+Sk )  the  names  2 x(n)  and  (n)  respectively,  with  notational 
dependence  on  time  n  restored  in  the  names.  The  name  is  not  meant  to  suggest  time  but  respects 
a  different  naming  convention.  The  factor  of  2  in  2 x(n)  gives  unity  DC  gain  to  the  filter  that  has 
x(n)  as  its  impulse  response,  and  it  also  puts  the  frequency-domain  version  of  spreading  equation 
([2])  into  a  standard  form.  In  the  time  and  frequency  domains  these  definitions  and  Q  become 


x(n)  4  5(n+l)  +  5(n  — 1) 
2 


Pj27r/i  p-j27r/ 

X{f)  = - i - =  cos(2tt/) 


(3) 


to(n)  =  8(n)  <-> 

t\{n)  =  x(n)  <-> 

4(n)  4  2x(n)  *4_i(n)  -  tk-2(n ) 


=  1 

T1(X(/))  =  X(/) 

Tk{X{f))  =  2X(f)Tk.1(X(f))  -  Tfc_2(X(/)),  (4) 


with  the  last  line  applying  for  k  >  1.  On  the  right  Tk(f)  has  been  written  Tk(X(f))  to  emphasize 
that  in  fact  it  depends  on  /  only  through  the  value  attained  by  X(f).  Further,  written  in  this  way 
each  Tk  is  in  fact  a  polynomial  in  its  argument.  When  the  /  dependence  of  X  (/  )  is  suppressed,  the 
first  few  polynomials  are 


To  =  1 

Ti=X 

T2  =  2X2  -  1 

T3  =  4X3  -  3X 

T4  =  8X4  -  8X2  +  1. 
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Fig.  3  Replacing  the 
default  one-add  X0\ d(/) 
with  a  five-add  Xn&N(f) 
(bottom,  right  scale) 
here  transforms  a  half¬ 
band  H(X0\d(f))  with 
4  nonredundant  coeffi¬ 
cients  into  H(Xnmi(f)) 
(top,  left  outside  scale, 
dB).  An  alternative 
X0\d  (/)  requires 
seven  more  coefficients 
than 

0.635 


A'0id(/) 

to 

perform 

simi- 

larly.  The 

plot 

(  inner  scale 

of 


(/old)  — 
old  (x  new  (/old)) 


/new 

^r1 

characterizes  the 
-  -0.635  net  warping  of 
frequency. 


0.1405  0.212  0.25  0.288 


0.3595 


These  polynomials  Tq,Ti,T2,  . . .  are  actually  the  Chebyshev  polynomials  (of  the  first  kind),  for 
which  the  name  T \  is  in  fact  customary.  Impulse  response  ([I])  can  be  written  in  both  domains  using 
these  Chebyshev  polynomials  and  their  time-domain  equivalents  as 

M  M 

h(n)  =  h0t0(n)  +  J22hktk(n)  <-+  =  h0  +  £  2hkTk(X(f)).  (5) 

k= 1  k=  1 


This  relationship  and  the  Chebyshev  recursion  above  can  be  realized  jointly  to  create  a  filter  in 
what  we  might  reasonably  call  Chebyshev  direct  form  as  shown  in  Fig.  [I]  The  reversed  Chebyshev 
form  (8)  of  Fig.  [2]  is  just  the  transpose  of  the  Fig.  [I]  filter,  the  latter  with  every  arrow  of  the  signal 
flow  graph  reversed.  The  component  impulse  response  x(n)  in  these  figures  is  as  defined  in  ([3])  and 


so  is  not  causal.  Realization  therefore  requires  shimming  delays  where  shown 


2.2  A  transformation:  replace  the  linear-phase  spreading  function  x 

The  temporal  structure  of  the  impulse  response  on  the  left  side  of  ([5])  is  completely  contained  in 
the  functions  to,  . . . ,  which  contain  impulses  at  various  times,  but  it  is  the  mapping  in  ([5])  from 
these  to?  >  •  •  •  to  h(n)  that  contains  all  the  coefficient  information  and  hence  the  essence  of  this 
filter’s  particular  frequency  response.  It  is  this  separation  that  is  interesting  and  useful,  and  it  takes 
its  cleanest  form  in  the  frequency  domain.  In  fact,  on  the  right  in  ([5])  the  frequency  response  is 
written  as  H(X (/))  rather  than  the  more  conventional  H(f)  to  emphasize  that  it  depends  on  /  only 
through  its  dependence  on  X(f).  In  the  two-step  mapping  f  —>  X  —>  H  in  effect  /  determines  the 
value  of  X (/),  which  then  determines  the  value  of  H(X (/)).  This  means  we  can  change  H(X (/)) 
not  only  in  the  usual  way,  by  changing  the  coefficients  h that  determine  the  second  step  X  — >  H, 
but  by  replacing  the  function  X(f)  with  a  different  function  so  as  to  affect  the  first  step  /  — ►  X. 
We  can  replace  it  not  only  in  analysis  but  in  the  implementation  itself. 
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Fig.  4  The  unit  cells  of  Fig.  [2]  (up¬ 
per  corners  here)  can  be  realized  as 
FPGA  unit  cells  (center)  realizing  the 
standard  linear-phase  spreading  func¬ 
tion  xoid  (n)  (top)  or  the  alternative 
^new(n)  of  Fig.  [§]  (bottom).  Circles 
are  adders,  rectangles  are  delays,  and 
block  delays  are  numbered.  Scalings 
by  two  and  by  1/2  are  A  and  V  respec¬ 
tively.  Colors  denote  special  cases:  for 
the  left  unit  cell,  for  the  right  unit  cell, 


Consider  the  example  in  Fig.  [3]  There  is  a  lot  of  information  here,  so  it’s  best  walked  through 
slowly.  Three  distinct  vertical  scales  are  used,  here  termed  the  left,  right,  and  inner  scales,  but  the 
horizontal  scale  is  the  same  for  all  curves  and  is  just  ordinary  normalized  frequency  /.  On  the  top 
the  magnitude  response  of  an  example  prototype  filter  H(X0\d(f))  of  length  15  is  plotted  against 
/  using  the  left  vertical  dB  scale.  This  particular  prototype  has  been  designed  as  a  halfband  filter, 
so  on  a  linear  scale  the  passband  and  stopband  ripples  are  mirror  images  and  the  response  is  down 
by  1/2  or  about  6  dB  at  /  =  1/4.  On  the  bottom  is  plotted,  against  /  on  the  linear  right  scale, 
X0id (/)  =  cos(2tt/),  the  Fourier  transform  of  the  default  linear-phase  spreading  function 

£oid  (n)  =  \  (<5(n+l)  +  5{n- 1)). 


Replacing  X0\&  (/)  with  Xnew(/),  the  Fourier  transform  of  modified  linear-phase  spreading  function 

^new(^)  =  ^(S(n+1)  +  S(n- 1)) 

-  |(5(n+3)  +  S(n—3)) 

+  ^(5(n+5)  +  8{n— 5)), 


replaces  prototype  response  H(X0\&(f))  with  transformed  response  H(Xnew(f))  in  these  plots 
by,  in  effect,  warping  the  frequency  axis.  It’s  as  if  at  any  particular  frequency  f0 id  the  values  of 
X0\d(foid)  and  #(^oid(/oid))>  one  above  the  other  in  the  plots,  are  linked  as  shown  by  the  right¬ 
most  of  the  two  vertical  lines  in  Fig.  [3]  for  an  example  /0id  •  The  link  is  then  moved  horizontally 
(left  in  this  example)  to  frequency  /new  such  that  X0id(f0\d)  =  Xnew(/new).  Then  H(X0id(fM))  = 

-ff  (^new(/new))  aS  Well. 

There  is  another  way  to  characterize  the  frequency  warping  graphically,  one  that  we  will  prefer 
in  subsequent  discussion.  Since  X0id (/)  is  strictly  decreasing  on  the  closed  interval  [0,0.5],  it  is 
invertible  on  that  interval.  Here  the  range  of  Xnew(/)  is  contained  in  (equals,  in  fact)  that  interval, 
so  V, id  (./old)  =  -View  (/new)  implies  /old  =  V/]  (View  (/new))  or,  for  this  XoU(f)  in  particular, 
/old  =  ^  acos(Xnew(/new))-  In  Fig.  lathis  /oW  is  plotted  faintly  as  a  function  of  /new  on  a  third, 
vertical  scale  and  with  a  •  wherever  this  plotted  /0id  crosses  a  horizontal  grid  line.  The  shape 


inner 


of  this  curve  reveals  the  character  of  the  warping.  A  line  with  unit  slope  would  represent  no  warping 
at  all.  Where  this  curve  is  steeper  or  shallower  than  that,  the  /new  frequencies  are  respectively  being 
squeezed  closer  together  or  stretched  further  apart. 

For  this  particular  filter-design  example,  a  preliminary  resource  count  for  Virtex-family  FPGA 
implementation  follows  easily  from  the  Fig.  [4]  sketches  of  pipelined  unit-cell  realizations  for  the  re¬ 
versed  Chebyshev  form  (which  may  or  may  not  be  more  efficient  than  the  Chebyshev  direct  form). 
In  the  Virtex  family  one  slice  comprising  two  lookup  tables  is  required  for  two  wordwidth  bits  of 
a  delayed  adder,  a  pure  delay,  or  a  block  delay.  Taking  zero  coefficients  into  account,  the  new 
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transformed  version  of  the  Fig.  [3]  halfband  filter  requires  four  nontrivial  coefficient  multiplies  and 
87  slices  per  wordwidth  bit  pair.  Obtaining  similar  performance  by  simply  redesigning  the  origi¬ 
nal  filter  coefficients  to  replace  prototype  response  with  the  roughly  similar  response 

X0id (/)  that  is  for  comparison  in  Fig.  [3] requires  1 1  nontrivial  coefficient  multiplies 

and  94  slices  per  wordwidth  bit  pair  when  realized  in  conventional  linear-phase  transposed  form. 
Relative  to  this  conventional  approach,  the  direct  realization  of  transformed  filter  has 

saved  seven  (nearly  2/3)  of  the  coefficient  multiplies  and  seven  slices  (about  7%)  per  wordwidth  bit 
pair  (and  it  achieves  slightly  better  performance). 

3  The  idea  applied  to  2D  filters:  the  McClellan  transformation 

While  these  sorts  of  efficiency-improvement  transformations  can  be  useful  in  ID  filter  design,  it’s 
in  2D  filter  design  that  the  technique  is  the  most  powerful.  In  this  context  it  is  termed  the  McClellan 
transformation  J2  •  Here  we  use  it  to  design  weights  for  shaping  the  beam  of  a  planar  antenna  array. 

There  are  several  stages  to  our  transition  to  the  2D  planar-array  application.  First,  we  restate 
our  earlier  discrete-time  mathematics  in  continuous-time  form.  Second,  we  replace  time  with  space 
and  move  from  ID  to  2D.  Finally,  we  explore  the  planar-array  application  through  examples. 

3.1  Replace  discrete  time  with  continuous  time 

Every  discrete- time  function  has  a  continuous-time  equivalent.  Using  a  bar  to  indicate  a  function  of 
continuous  time,  the  continuous-time  equivalent  of  x(n),  for  example,  is 

m  =  E  x(n)S(t  —  Sn).  (6) 

n 

The  sample  spacing  would  more  typically  be  called  T,  but  that  letter  has  quite  enough  uses  here 
already,  so  S  for  “Spacing”  has  been  chosen  instead. 

If  we  likewise  use  the  bar  to  denote  the  continuous-time  frequency  variable  /,  the  continuous¬ 
time  Fourier  transform  of  ([6])  is 

f  (j2x(n)  5(t  -  Sn)\ e~i2nIt  dt  =  J2  x(n )  e“j27r/Sn  =  X(fS), 

^  '  n  '  n 

where  X  (/)  is  the  discrete-time  Fourier  transform  used  before.  We  get  the  continuous-time  trans¬ 
form  simply  by  evaluating  discrete-time  transform  at  /  =  fS.  We  will  write  all  continuous-time 
transforms  this  way  and  so  will  speak  of,  for  example,  continuous-time  Fourier  pair  x(t)  X  ( fS ). 

Of  course  then  the  product  X(fS)Y(fS )  of  continuous-time  transforms  is  just  the  product 
X(f)Y(f )  of  discrete-time  transforms  evaluated  at  f  —  fS,  implying  that  x(t)*y(t),  which 
uses  a  continuous-time  convolution,  is  the  continuous-time  equivalent  of  x(n)  *y(n ),  which  uses  a 
discrete-time  convolution.  Equivalence  ([6])  is  also  linear  in  x(n),  so  any  time-domain  relationship 
involving  discrete-time  functions,  scaling  by  constants,  addition  and  subtraction  of  functions,  and 
convolution  of  functions  has  a  continuous-time  equivalent  that  can  be  obtained  by  simply  replacing 
its  constituent  functions  with  continuous -time  equivalents  individually  and  then  taking  all  convolu¬ 
tions  to  be  in  continuous  time.  It  follows  that  the  continuous-time  equivalents  of  the  definitions  and 
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relationships  of  0  and  ([5])  are  just 

t0(t)  =  6  (t )  <-> 

Tq(X)  =  1 

t\(t)  =  x{t)  <-> 

Tl(X)  =  X 

tk(t)  =  2x(t)*tk_1(t)  -  tk-2(t )  ^ 

Tk(X)  =  2XTk_1(X)-Tk_2(X), 

M 

M 

h(t)  =  h0t0(t)  +  y^2  hktk(t) 

H(X(fS ))  =  h0  +  2  hkTk(X(fS)). 

k= 1  k= 1 


Each  occurrence  of  X(fS ))  in  the  recursion  on  the  right  has  been  written  as  just  X  for  brevity 
and  to  emphasize  that  the  form  of  that  recursion  is  precisely  that  of  the  classic  description  of  the 
Chebyshev  polynomials. 

The  important  thing  here  is  that  we  still  have  a  two-step  frequency-domain  mapping,  now  /  — ► 
X  — ►  H.  We  can  modify  H(X(fS ))  by  changing  filter  coefficients  to  change  the  second  step 
X  — >  H  or  by  replacing  the  linear-phase  spreading  function  x(t)  to  change  the  first  step  f  —>  X. 


3.2  The  planar-array  application 


The  2D  FIR  impulse  and  frequency  responses  of  most  interest  to  us  are,  respectively,  array  tapers 
and  array  factors,  and  here  we  define  them  in  a  way  that  makes  (array  taper)  (array  factor)  a 
Fourier  pair  in  various  ways  corresponding  to  different  choices  of  normalization. 

The  usual  periodic  structure  of  the  antenna  elements  of  a  narrowband  planar  antenna  array  can 
be  characterized  by  nominal  element  locations  that  are  integer- weighted  combinations  of  3D  basis 
vectors,  call  them  Asi  and  As2,  where  normalizing  wavelength  A  =  c/|/|  for  |/|  at  the  lower  edge  of 
the  signal  band.  If  dimensionless  vectors  s±  and  S2  are  made  the  columns  of  3  x  2  element- Spacing 
basis  matrix  S,  every  element  then  has  a  nominal  location  of  form  ASn,  where  element  index  n  is 
a  two-element  integer  column  vector  that  identifies  the  element. 

We  can  arrive  at  a  Fourier  (array  taper)  (array  factor)  relationship  by  associating  element 
weights  to  element  positions  normalized  in  one  of  four  ways: 


weight  h{n)  associated  with  position  < 


n 

S  n 

±XS  n 


fully  normalized, 
normalized  by  A, 
no  normalization. 


In  the  formulations  that  are  most  convenient  we  take  h(n)  to  be  the  weight  associated  with  the 
element  at  physical  location  ASn,  but  we  can  choose  whether  to  carry  full  position  information 
throughout  or  to  instead  let  either  A  or  AS  remain  unstated  side  information. 

Fully  normalized  array  taper  h(n)  can  be  taken  to  be  a  2D  discrete-“time”  (here  actually  space) 
function  with  2D  discrete-“time”  (space)  Fourier  transform 

H(f)  =  J2h(n)  e-^n  .  (8) 

n 

This  sum  and  others  like  it  to  follow  are  over  all  possible  2D  element  indices  n,  but  of  course 
element  weight  h(n)  is  nonzero  only  for  those  n  corresponding  to  actively  used  elements.  Here 
(array  taper)  ^  (array  factor)  takes  the  form  h{n)  ^  H(f ),  a  2D  Fourier  pair  that  makes  the 
unit  square  in  2D  normalized  frequency  f  a  period  of  fully  normalized  array  factor  H(f).  Of 
course  Fourier  pair  h{n)  H(f)  mirrors  exactly  the  structure  of  a  2D  FIR  filter  and  so  motivates 
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this  entire  report,  which  is  fundamentally  about  applying  the  McClellan  transformation,  a  classic 
technique  for  2D  FIR  filter  design,  to  the  design  of  array  factors. 

An  array  taper  normalized  by  A  takes  the  form 

My)  =  MM  My  -  SrM  (9) 

n 

where  3D  column  vector  Ay  of  reals  is  position  and  y  itself  is  therefore  position  normalized  by  A. 
Element  positions  are  represented  by  impulses  with  element  weights  as  areas.  Substituting  ([8])  and 
([9])  on  the  left  and  right  readily  verifies  the  equality  in 

H(rS)  =  J  h(y)  e~^ay  dy,  (10) 

where  the  integral  is  in  3D  with  dy  a  volume  differential  and  where  a  is  a  real-valued  3D  row 
vector.  This  provides  an  (array  taper)  ^  (array  factor)  relationship  in  the  form  h( y)  ^  iT(aS),  a 
continuous-“time”  (space)  3D  Fourier  pair  with  spatial  frequency  a  as  the  transform  variable.  The 
curious  fact  that  3D  function  H( aS)  can  be  constructed  from  2D  function  H(f)  is  fundamentally 
due  to  the  array  being  confined  to  a  plane,  a  2D  subset  of  3D. 

The  array  factor  in  this  A-normalized  form  has  natural  interpretations.  For  a  receive  array,  the 
component  of  the  array  output  arising  from  signals  with  direction  of  arrival  (DOA)  given  by  3D 
row  unit  vector  a  is  just  H(a S)  times  the  unweighted  output  of  the  single  element  at  the  spatial 
origin^  For  a  transmit  array  it  is  a  little  more  complicated  @,  but  if  we  assume  that  the  far- field 
radiated  electromagnetic  field  due  to  a  signal  input  to  an  element  is  everywhere  linear  in  the  product 
of  that  signal  with  the  weight  associated  with  that  element,  the  electric  or  magnetic  field  at  a  point 
in  the  far  field  in  direction  a  is  then  just  H(a S)  times  the  value  of  that  same  field  when  a  single 
unweighted  element  at  the  spatial  origin  is  used  alone  to  transmit. 

There’s  actually  more  than  one  way  to  state  the  (array  taper)  (array  factor)  relationship  with 
“no”  normalization,  because  a  free-space  plane  wave  can  be  expressed  with  equal  ease  in  terms  of 
a  radian-frequency  wavenumber  vector,  ^  times  that  same  vector,  or  the  negative  of  either.  For 
receive  arrays  that  negative  vector  points  in  the  DOA  and  so  is  quite  natural.  Most  natural  for  us 
is  to  use  spatial-frequency  row  vector  k  =  which  makes  the  transmit-array  and  receive-array 
radian  wavenumber  vectors  respectively  27rk  and  —  27rk.  Using  unnormalized  position  x  =  Ay  and 
spatial  frequency  k  =  then,  ([8])  yields 

ff(kAS)  =  MM  e“j2^kASn  =  f  (j2h(n)  5 (x  -  ASn)  j  e“j27rkx  dx,  (1 1) 

n  **  '  n  ' 

thereby  establishing  Fourier-pair  relationship 

J2h(n)6(x-  ASn)  <-»  H{ kAS) 

n 

between  variables  x  and  k.  (As  a  derivation  strategy  this  direct  integration  of  the  delta  functions  in 
(|TT|)  is  simpler  than  carrying  out  a  change  of  variables  on  the  integral  in  (fTO]),  because  scaling  the 
independent  variable  in  a  delta  function  is  a  little  nonintuitive  and  for  most  people  requires  special 
care.) 

As  an  aside,  note  that  technically  k  can  just  as  well  be  negated  relative  to  this  so  that  radian 
wavenumbers  for  a  receive  array  take  form  27rk,  but  that  in  this  case  the  arguments  to  both  h 

'The  length  of  this  “unit  vector”  is  precisely  one  only  at  the  bottom  of  the  signal  band.  For  a  truly  narrowband  array 
of  course,  that  is  all  that  matters. 
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and  H  become  negated  as  well.  This  makes  h(—n )  the  weight  associated  with  element  position 
AS n  and  thereby  moves  the  awkward  minus  sign  from  the  spatial-frequency  domain,  where  it  only 
meant  thinking  of  DOA  instead  of  propagation  direction,  to  the  time  domain,  where  there  is  just  no 
convenient  interpretation.  This  minus  sign  can  be  moved  to  many  places  in  the  formulation,  but  in 
the  end  it  will  not  go  away. 

Overall  then,  the  “fully  normalized”  formulation  best  ties  the  array  application  to  FIR  filter 
design,  but  the  “normalized  by  A”  formulation  with  position  =  Sn  will  be  the  default  going  forward 
because  unlike  the  fully  normalized  formulation  it  retains  array  size  and  shape  information,  makes 
classic  element  spacings  like  A/2  and  A/ a/3  easy  to  see,  makes  common  array-taper  symmetries 
readily  visible,  and  maps  azimuth  and  elevation  coordinates  onto  array-factor  plots  in  an  intuitive 
way. 

3.2.1  Steering  means  we  cannot  restrict  to  unit  vectors 

Discrete-“time”  Fourier-transform  property  h(n)  e^27rsSn  ^  H(f  —  s' S)  in  2D  implies  that  using 
weights  phase  shifted  in  this  way  in  the  array  taper  results  in  an  array  factor  H^a  —  's) S)  in  which 
features  near  the  origin  in  H(a S)  appear  instead  near  some  DOA  steering  vector  s'.  Given  this 
background  fact,  it  is  not  necessary  to  explicitly  consider  steering  here,  except  to  realize  that  while 
a  and  s'  are  each  unit  vectors,  their  difference  a  —  s'  in  general  is  not.  That  difference  vector  may 
have  a  length  approaching  two.  In  respect  of  that,  we  generally  compute  and  plot  iT(aS)  against 
3D  row  vector  a  without  requiring  vector  a  to  be  a  unit  vector.  Of  course  any  boresight  component 
in  vector  a  is  orthogonal  to  both  basis  vectors  and  therefore  has  no  effect  on  product  aS,  so  it  is 
enough  to  describe  the  behavior  of  H( aS)  versus  the  array-plane  component  of  vector  a. 

3.2.2  What  array  tapers  are  desired? 

We  will  always  design  our  unsteered  iT(aS)  to  have  a  main  beam  at  boresight,  which  is  just  the 
line  of  a  values  orthogonal  to  the  array  plane  and  therefore  to  both  si  and  S2  so  that  f  =  aS  =  0. 
This  means  we  want  iT(aS)  to  be  a  lowpass  filter  in  the  array-plane  component  of  spatial  frequency 
a.  Most  often  we  want  a  narrow  passband,  and  quite  commonly  we  are  willing  to  allow  iT(aS) 
to  display  various  kinds  of  symmetry,  though  we  may  not  actually  need  such  symmetry.  Typical 
symmetries  of  this  type  are  rotations  about  boresight  and  reflections  through  lines  in  the  array  plane 
that  pass  through  the  origin. 

Array  tapers  typically  have  linear  phase,  yet  another  form  of  symmetry,  though  it  is  in  no  way 
necessary  in  general.  It  is,  however,  necessary  if  the  McClellan  transformation  explored  below  is 
to  apply.  Typically  the  spatial  origin  is  placed  at  the  point  of  symmetry  so  that  the  (array  taper)  ^ 
(array  factor)  in  fact  has  zero  phase:  h(y)  —  h*(— y)  or  h(n)  =  h*(—n)  with  H( aS)  real. 

Unsteered  element  weights  are  often  real  valued  as  well  so  that  both  the  unsteered  array  taper 
and  the  unsteered  array  factor  are  real  and  even.  This  will  be  true  of  all  the  design  examples  in  this 
report,  even  though  the  McClellan  transformation  does  not  in  any  way  require  it. 

3.2.3  Taper  loss 

One  simple  characterization  of  a  receive  array  is  boresight  SNR,  the  ratio  of  the  array  output  power 
due  to  a  plane-wave  signal  arriving  from  the  boresight  direction  to  the  noise  power  output  by  the 
array.  Suppose  all  array  elements  have  identical  but  separate  preamps  with  enough  gain  that  their 
internally-generated  noise  levels  dominate  array  output  noise.  In  this  case  these  array-output  signal 
and  noise  power  levels  are  respectively  proportional  to  |i7(0)|2  and  to  \\h\\ 2  =  \hn\2,  the  sum 

over  all  elements  of  the  squared  magnitudes  of  the  weights. 
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One  simple  characterization  of  a  transmit  array  is  the  ratio  of  the  far-held  signal  power  per 
steradian  radiated  in  the  boresight  direction  to  the  total  power  delivered  to  the  elements  by  the 
transmitter.  Suppose  the  transmitter(s)  is(are)  well  matched  to  the  elements  at  all  steering  angles.  In 
this  case  this  radiated  signal-power  density  and  the  power  delivered  to  the  array  by  the  transmitter 
are  proportional  to  |  H (0)  | 2  and  to  1 1  h  \  | 2  respectively. 

So  for  both  transmit  and  receive  arrays  the  ratio  \H(0)\2  /\\h\\2  is  meaningful  as  a  performance 
measure.  We  don’t  know  the  constants  of  proportionality  that  would  enable  us  to  directly  interpret 
these  numbers,  but  we  can  sidestep  this  issue  by  comparing  this  measure  of  array  performance  to 
the  performance  of  a  reference  array  described  by  discrete-time  Fourier  pair  r(n)  R(f).  This 
report  follows  common  practice  and  takes  this  reference  array  to  be  N  uniformly  weighted  elements 
so  that  r(n)  =  1  and  R( 0)  =  ||r||2  =  N. 

In  the  N  =  1  special  case  corresponding  to  a  single  embedded  element  used  alone  as  a  one- 
element  array,  the  performance  ratio  R(0)/\\r\\2  =  N  goes  to  unity.  In  the  more  general  case 
this  ratio  i?(0)/||r||2,  for  the  reference  array,  or  H(0)/\\h\\2,  for  any  array,  can  be  interpreted  as  a 
performance  gain  relative  to  this  single-element  case.  For  a  receive  array  this  is  classic  SNR  gain , 
and  for  a  transmit  array  it  is  classic  power  gain. 

While  these  gain  figures  are  often  given  explicitly,  in  this  report  we  will  report  them  only  relative 
to  reference-array  gain: 

|ff(0)|2/IHI2  |ff(0)|2/IHI2  lEn^l2  \(h1r)l 

|i?(0)|2/||r||2  N  \\h\\2  N  ||/i||2||r||2, 

using  on  the  right  inner  product  |  (/i,  r)  |2  =  and  squared  norm  \\h\\ 2  —  (h,h).  The  point  of 

this  seemingly  odd  formulation  is  that  by  the  famous  mathematical  inequality  of  Cauchy,  Schwarz, 
and  Bunyakovsky  in  Hilbert  spaces  the  ratio  on  the  right  cannot  exceed  unity.  This  implies  that  no 
taper  can  improve  on  the  gain  of  a  uniformly  weighted  array.  For  each  design  then  the  dB  distance 
below  this  gain  bound  will  be  presented  as 

I  (h  v)  1 2 

(taper  loss)  =  -101og10  p^p^p- 

While  one  generally  desires  a  small  taper  loss,  other  measures  of  performance  are  important  as 
well.  For  example,  one  might  want  the  best  sidelobe  rejection  possible  (by  any  of  several  measures) 
consistent  with  some  tolerable  level  of  taper  loss. 


3.2.4  Simple  basis  matrices  for  radar  arrays 

Here  are  clean  ways  to  get  the  basis  matrices  of  the  most  interest  in  radar,  assuming  we  can  orient 
the  array  plane  relative  to  the  coordinate  frame  in  any  way  we  like. 


Square  grid:  Choose  orthogonal  vectors  that  are  each  zero  in  two  components  and  of  unit  mag¬ 
nitude  in  the  third.  Scale  the  resulting  matrix  by  1/2.  Twenty-four  distinct  matrices  of  this 
type  are  possible,  and  they  describe  three  distinct  element  grids  with  the  optimal  A/2  spacing. 


Example:  \ 


o  o 
0  -1 
Ll  0J 


(This  one  has  array  boresight  in  the 


l 

o 

Lo  J 


direction.) 


Alternatively,  zero  any  one  common  component  of  two  noncollinear  vectors  and  give  the 
other  components  of  each  unit  magnitude.  Scale  the  resulting  matrix  by  1/ (2\/2).  This  yields 
another  24  distinct  matrices  describing  three  more  distinct  grids  with  A/2  spacing.  Example: 


l 

2V2 


1  -1 
0  0 
Ll  l 


(This  one  has  array  boresight  in  the 


direction.) 
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Equilateral  triangular  grid:  Give  the  two  vectors  a  single  zero  component  each,  and  place  that 
zero  in  different  components  for  the  two  vectors.  For  each  vector  give  each  of  the  two  other 
components  unit  magnitude.  Scale  the  resulting  matrix  by  1/V6.  This  yields  96  distinct 
matrices  describing  four  distinct  grids  with  the  optimal  A/\/3  spacing.  Half  the  matrices 


have  vectors  60°  apart  and  half  have  them  120°  apart.  Example: 


vectors  120°  apart  and  has  array  boresight  in  the 


l 

-l 

l 


l  -l 

1  o 
o  1 


(This  one  has 


direction.) 


3.3  The  transformation  in  2D 

The  Chebyshev-form  FIR  filter  structure  in  ID  continuous  time  with  impulsive  samples  spaced  by  S 
is  described  by  Q  and  (|7]>  and  has  frequency-response  structure  /  — >  Xm(fS)  — ►  H(Xm(fS)), 
with  the  latter  depending  on  scalar  continuous-time  frequency  /.  For  the  planar-array  application 
we  use  alternative  structure  — ►  X2d(  )  — ►  H(X 2d(  ))  in  which  the  second  constituent  map 

X  —>  H  is  exactly  as  before  but  is  now  preceded  by  the  mapping  — >  X2d(  )  of  a  3D  row- vector 
frequency  into  a  scalar- valued  X2d(  )  using  2D  FIR-filter  frequency  response  X^d  (/)  in  dis¬ 

crete  “time”  (space)  and  3x2  element- spacing  basis  matrix  .  The  original  ID  linear-phase  spread¬ 
ing  function  represented  by  continuous-time  Fourier  pair  xyo (t)  Xm(fS)  =  cos(2 nfS)  is 

replaced  in  the  planar- array  application  with  a  continuous-“time”  (space)  transform  pair  x^(  ) 
X2d(  ),  where  the  latter  3D  Fourier-pair  relationship 

^2d(  )  =  J  xM  )  e-j2T  d  (12) 

follows  from  definitions  that  determine  both  members  of  the  pair  from  coefficients  £2d(^)  indexed 
by  a  two- vector  n  of  integers: 


*2D(/)  ^  D(n)  e~'^n 

Tl 

(13) 

x3D(  )  =  y^X2D(«)<5(  -  n). 

(14) 

n 


Here  (p"3])  makes  X2d(ti)  ^  ^2d(/)  a  2D  discrete-time  Fourier  pair.  Relationships  (fl2|)  and  (fl4|), 
which  move  that  pair  relationship  into  3D  spatial  coordinates,  are  structurally  analogous  to  ([l0|) 
and  ([9]),  so  technically  3D  Fourier  pair  xi,d{  )  ^  X2d(  )  has  an  (array  taper)  <-+  (array  factor) 
structure  exactly  as  if  X2D (n)  were  the  element  weights  for  a  (presumably)  tiny  planar  array.  But 
of  course  there  is  no  such  array.  Instead,  the  utility  of  xi,b(  )  ^2d(  )  is  in  permitting  the  easy 

construction  of  a  (array  taper)  ^  (array  factor)  pair  for  a  large  array. 

That  construction  of  a  (array  taper)  <-+  (array  factor)  pair  was  effected  originally  through  the 
Chebyshev  recursion  in  ([7]),  which  used  repeated  time-domain  convolution  with  xyo(t)  to  spread  a 
seed  function  to  (t)  —  S (t)  farther  and  farther  out  in  time  to  provide  basis  functions  to  (t ) ,  t\ (t) ,  1 2  (t) , . 
from  which  an  arbitrary  linear-phase  filter  impulse  response  could  be  built.  The  new  Chebyshev  re¬ 
cursion  uses  repeated  convolution  in  the  array-plane  subspace  of  3D  “time”  (space)  in  a  similar  way. 
Convolution  with  ^3d(  )  spreads  a  seed  function  to(  )  =  S(  )  farther  and  farther  out  in  the  array 
plane  to  provide  basis  functions  to  (  ) ,  i\  (  ) ,  t*2  (  ) , . . .  from  which  an  arbitrary  linear-phase  array 
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taper  can  be  built.  The  planar-array  analog  of  0, 

i0(  )  =  S(  )  <-► 

/ 1 (  )  -  ®3d(  )  ^ 

h{  )  =  2. x3D(  )*4-i(  )  -  4-2 (  )  <-> 

M 

h(  )  =  Mo(  )  +  D2Mfc(  )  «-► 

fc=l 


T0(X)  =  1 
Ti(X)  =  X 

Tk(X)  =  2XTk_1(X)-Tk_2(X), 

M 

H(x 2D(  ))  =  /i0  +  ^2/lfcTfc(X2D(  )),  (15) 

/c=l 


follows  from  the  new  definitions  of  t o  (  )  and  X3D  (  ) ,  with  only  the  latter  representing  a  meaningful 
change  as  to(  )  is  and  must  be  simply  the  identity  element  for  convolution  in  the  function  space  in 
which  ^3d(  )  is  defined.  The  rest  is  derived  from  these  two  quantities  just  as  before. 

Because  we  need  X2d{  )  to  be  real  valued,  ^3d(  )  should  be  a  zero-phase  impulse  response, 
which  is  to  say  that  it  must  have  conjugate  symmetry  x 3d  (  )  =  £  3d  (—  )  or,  equivalently,  X2D (n)  = 
x 2d ( — xt).  We  also  need  the  range  of  X2d(  )  to  be  contained  in  the  closed  interval  [—1, 1],  which 
was  the  range  of  the  original  Xm  (/)  in  ID,  but  this  is  easy  to  arrange:  we  can  scale  and  offset  any 
preliminary  choice  of  a  zero-phase  X2d(  )  by  constants  as  necessary  to  fix  its  range  to  [—1, 1].  This 
scales  X3D  (  )  and  offsets  it  by  a  delta  function. 


3.3.1  Design  example:  the  simplest  spreading  function  has  seven  samples 

Let  us  work  with  a  concrete  example  as  we  continue.  For  now  suppose  in  particular  that  the  basis- 
vector  columns  of  matrix  have  length  1  /  y/3  and  that  they  are  60°  apart  (whether  or  not  chosen  in 
the  way  described  earlier  in  Section|T2.4|).  This  puts  the  array  elements  on  an  equilateral-triangular 
grid,  a  scaled  version  of  the  so-called  hexagonal  lattice,  with  elements  spaced  at  the  classic  maxi¬ 
mum  if  grating  lobes  are  to  be  avoided  at  all  scan  angles.  Array  taper  h(  )  will  comprise  impulses 
on  points  of  the  lattice  as  desired  if  and  only  if  ^3d(  )  comprises  such  impulses.  To  the  extent  that 
we  can  do  so  while  staying  in  the  lattice,  let  us  give  X3d(  )  roughly  equal  extent  in  all  directions. 
The  simplest  initial  choice  locates  unit  impulses  at  the  points  of  a  hexagon  centered  on  the  origin, 
at  the  origin’s  nearest  neighbors  in  the  lattice: 

®3d(  )  = 

—  [oD  +  <5(  —  [?])  +  5(  —  [-}]) 

+  —  [~o])  +  <5(  —  [_?])  +  <)(  —  [_!])• 

This  Fourier  transforms  to  an  array  factor  X2d(  )  comprising  three  cosine  terms.  It  has  range 
[—3, 6],  so  scaling  by  2/9,  which  yields  range  [  —  |] ,  and  subsequent  addition  of  constant  —1/3 

fixes  the  range  at  the  desired  [—1,1].  The  final  form  is 

®3D(  )  = 

-H  -  [oD 

+  I'K  —  [oD  +  §<H  —  [iD  +  §<5(  —  [~i]) 

+  §<H  —  [-oD  +  §£(  —  [-?])  + 1^(  —  [-!])•  (16) 

The  first  seven  Chebyshev-derived  basis  functions  in  the  “time”  (space)  domain  that  result  from 
this  -On  I  )  are  depicted  in  Fig.  [5] in  seven  separate  sketches.  Each  sketch  lies  entirely  in  the  array 
plane,  which  has  axes  shown  faintly.  The  top  sketch  represents  the  single  unit  impulse  of  fo(  )• 
Moving  down,  subsequent  sketches  depict  / 1  (  ) . t,,(  ) .  The  numerator  of  each  integer-ratio 
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Fig.  5  Impulse  _i_ 
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impulse  area  is  shown  above  its  impulse-location  dot  •, 
and  the  denominator  common  to  the  dots  of  that  sketch 
appears  to  the  right  of  the  sketch.  Basis  function  £*.(  ) 
has  impulses  on  a  hexagonal  subset  of  the  lattice  with 
k+ 1  points  on  a  side  and  containing  l+3k(k+l)  points. 
The  latter  ties  array  size  to  prototype  filter  length  when 
the  McClellan  transformation  is  used  with  this  partic¬ 
ular  linear-phase  spreading  function  ^3d(  ). 

It  was  apparent  from  the  outset  that  the  hexago¬ 
nal  shape  of  £30  (  )  would  result  in  basis  functions 
t\ (  ) ,  £2  (  ) ,  •  •  •  with  support  confined  to  a  hexagonal 
subset  of  the  lattice,  and  indeed  the  nature  of  convo¬ 
lution  often  allows  one  to  anticipate  the  general  shape 
of  the  basis-function  support.  The  impulse  areas  are 
another  matter,  however.  The  pattern  of  impulse  ar¬ 
eas  in  Fig.  [5]  quickly  becomes  nonintuitive  as  k  in¬ 
creases.  The  )  sketch  at  the  bottom  even  shows  12 
areas  as  zero:  those  impulses  are  simply  absent.  The 
symmetries  in  the  areas,  the  invariance  to  various  ro¬ 
tations  and  reflections,  are  certainly  no  surprise  given 
the  presence  of  like  symmetries  in  x 30  (  ) ,  but  the  only 
viable  approach  to  learning  anything  more  about  the 
impulse  areas  appears  to  be  to  simply  use  the  recur¬ 
sion  on  the  left  side  of  ([15])  to  compute  them,  as  was 
done  here. 

Implementing  such  a  recursion  to  compute  either 
the  basis  functions  or  the  array  weights  themselves 
is  not  difficult.  Figure  [6]  presents  a  matlab  realiza¬ 
tion  of  a  2D  taper  computation  integrated  into  the 
basis-function  recursion  by  following  the  reversed- 
Chebyshev-form  FIR-filter  structure  of  Fig.  [2]  The 
Fig.  [6]  code  as  shown  operates  on  ID  Taylor  weights 
and  the  simple  x^(  )  of  ([16]),  but  any  ID  linear-phase 
taper  or  FIR  weights  can  be  substituted  for  the  former, 
and  any  linear-phase  x 3d  (  )  structured  as  impulses  on 
a  lattice  can  be  substituted  for  the  latter.  McClellan 
transformation  function  ftrans2  from  matlab’s  Im¬ 
age  Processing  Toolbox  should  provide  an  alternative 
to  most  of  the  Fig.  [6] code.  However,  this  author  is  not 
a  user  of  that  toolbox  and  so  has  not  experimented  with 
that  function. 

Figure  [8]  illustrates  the  transformation  of  an  exam¬ 
ple  ID  array  factor  into  2D  using  the  simple  X3d(  )  of 
([16])  in  a  McClellan  transformation.  The  array-plane 
contour  plot  on  the  left  of  the  second  row  shows  sev¬ 
eral  items  together.  The  plane  is  viewed  from  “behind” 
along  the  boresight  direction,  and  the  origin  is  in  the 
center. 
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w— 1 ; M=3  9 ; h=taylorwin (2  *M+1 , 12,-35) ; h=h/ sum (h)  ; 
x— [  0  2  2  ;  .  .  . 

2-3  2  ;  .  .  . 

2  2  0 ] / 9 ; cx=2 ; 

t=cell(l,3);  ct^zeros  ( 1 , 3 ) ; 
older=l;  t{older}=0;  ct (older) =1; 
old=2 ;t{old}=2*h(l) *w;  ct (old)=l; 
new=3 ; 
for  m=2:M 

tot=conv2 (2*x,t{old}) ; 
c=cx+ct (old) -1;  d=ct (older) -1; 
tot  (c, c) =tot (c, c) +  2*h (m) *w; 

tot (c-d : c+d, c-d : c+d) =tot (c-d : c+d, c-d : c+d) -t { older } ; 
t{new}=tot;  ct(new)=c; 

next=older;  older=old;  old=new;  new=next; 
end 

out=conv2 (x, t { old} ) ; 

c=cx+ct (old) -1;  d=ct (older) -1; 

out  (c,  c)  =out  (c,  c) +h  (M+l )  *w; 

out (c-d: c+d, c-d: c+d) =out (c-d: c+d, c-d: c+d) -t {older } ; 
out=out *  ( 1+3*M* (M+l ) ) ;  %optional  scaling  for  this  x 


Fig.  6  Matlab  realization  of  the  McClellan 
transformation  computes  w  h(  )  from 

•  scalar  input  w,  here  unity, 

•  ID  linear-phase  weight  vector  h  of 
length  2m+1,  here  the 

from  the  top  line  of  Fig. [8]  and 

•  2D  linear-phase  spreading  function 
^3d(  ),  here  from  example  {T6]i. 

The  reversed-Chebyshev-form  iteration  of 
Fig.  [2]  is  used  with  each  2D  function  rep¬ 
resented  by  a  square  matrix  and  the  index 
of  its  center.  For  example,  x3D(  )  here  be¬ 
comes  matrix  x  and  center  index  cx  with 
x ( i ,  j )  representing  x(i,j)(5(  -  [££]). 
The  code  assumes  m  >  1.  Result  w  h(  ) 
appears  as  matrix  out  and  center  index  c. 


Basis  vectors  are  shown  in  the  center  with  the  first  and  second  columns  of  basis  matrix  in  different 
colors. 

Linear-phase  spreading  function  )  is  shown  as  impulse-location  marks,  index  vectors  [^] 
that  specify  those  locations  in  ["*]  form,  and  the  impulse  areas  in  a  legend  at  the  lower 
left.  These  are  intentionally  made  a  bit  hard  to  see,  so  as  not  to  distract  from  the  main  items 
described  below. 

Having  the  first  basis  vector  point  downward  to  the  right  and  the  second  point  rightward  allows 
these  directions  to  correspond  to  how  row  and  column  indices  increase  on  the  page  when  a  weight 
matrix  is  presented  with  indentation  increasing  for  each  successive  line,  in  the  style  of  x  in  Fig.  [6] 
The  scale  of  the  two  plot  items  above  is  implied  by  the  earlier  discussion  of  the  basis  vectors 
used.  (Recall  that  they  have  been  given  length  l/y/3  here  to  indicate  element  spacing  of  A/\/3.) 
These  next  items,  in  that  same  plot  box,  are  on  a  completely  separate  and  unrelated  scale  from  the 
above.  They  share  the  box  with  the  above  only  to  save  space. 

A  projected  hemisphere  with  longitude  and  latitude  lines  at  10°  intervals  respectively  representing 
azimuth  and  elevation  relative  to  a  boresight  vector  normal  to  the  array  plane.  Here  the 
array  plane  and  the  boresight  vector  are  assumed  to  be  oriented  vertically  and  horizontally 
respectively,  and  azimuth  increases  in  the  clockwise  direction  when  the  array  is  viewed  from 
above.  This  hemisphere  has  unit  radius,  so  alternatively  one  can  view  it  as  a  3D  globe  with 
direction  unit  vector  a  extending  from  its  center  to  its  surface.  The  viewer  is  in  effect  in  the 
center  of  the  globe  looking  through  the  page  at  the  hemisphere  bulging  out  on  the  far  side. 

The  hexagonal  boundary  of  a  period  of  array  factor  X2d(  )•  It  circumscribes  the  hemisphere. 
Of  course  frequency- warping  function  X1d1(X2d(  ))  and  array  factor  H(X 20  (  ))  have 

the  same  periodicity  as  X2d(  )• 

The  frequency-warping  function  X1d1(X2d(  ))  =  ^  acos(X2D(  ))  in  contour-plot  form  as 

a  function  of  the  array-plane  component  of  vector  .  See  Table |T| for  contour  values. 

The  top  row  of  Fig.  [^presents  the  dB  magnitude  of  ID  prototype  response  H(Xm(fS))  with 
the  vertical  axis  extending  from  —60  dB  to  1  dB,  with  tick  marks  at  multiples  of  10  dB,  and  with  the 
horizontal  axis  representing  fS  values  from  0  to  0.5  if  the  response  is  interpreted  as  a  linear-phase 
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Table  1  The  contours  of  the  frequency-warping  plots  correspond  to  horizontal-axis  values  from  the  trans¬ 
formed  1 D  response.  That  axis  can  be  interpreted  as  the  usual  normalized  discrete-time  frequency  of  an  FIR 
filter  or  as  the  direction  cosine  of  a  classical  A/2  spaced  uniform  line  array. 


normalized 

FIR  frequency  fS 

direction  cosine  A / 
assuming  S'  =  A/2 

in 

plane 

major 

contours 

minor 

contours 

major 

contours 

minor 

contours 

at  the 

0 

0 

origin 

0.02  0.04  0.06  0.08 

0.04  0.08  0.12  0.16 

0.1 

0.2 

0.12  0.14  0.16  0.18 

0.24  0.28  0.32  0.36 

0.22  0.24  0.26  0.28 

0.44  0. 

48  0.52  0.56 

0.3 

0.6 

0.32  0.34  0.36  0.38 

0.64  0.68  0.72  0.76 

0.4 

0.8 

near 

peripheral 

extremes 

0.42  0.44  0.46  0.48 

0.84  0.88  0.92  0.96 

0.5 

1.0 

FIR-filter  response  and  representing  A /  values  from  0  to  1  if  the  response  is  interpreted  as  the  array 
factor  of  a  uniform  line  array  with  classic  element  spacing  S  =  A/2  and  a  linear-phase  taper.  In  the 
latter  case  axis  value  A /  is  the  direction  cosine,  the  cosine  of  the  angle  from  the  line  of  the  array  or 
the  sine  of  the  angle  from  boresight.  The  two  interpretations  are  offered  here  and  in  Table [ljbecause 
the  filter  and  line-array  views  will  likely  be  more  natural  to  signal-processing  and  array  researchers 
respectively.  Horizontal  axis  values  corresponding  to  contour  levels  are  marked  with  short  vertical 
lines  just  above  the  plot’s  centerline. 

The  top-row  prototype  response  is  actually  shown  for  two  choices  of  its  79  weights  h(n).  The 
,  wider  line  corresponds  to  Taylor  weights  with  a  —  35  dB  sidelobe  level  and  n  =  4.  Taylor 
weights  are  not  optimal  in  any  way  but  may  be  convenient  for  designers  lacking  optimization  tools 
or  the  time  to  deal  with  the  learning  curves  of  such  tools.  The  dark,  narrower  line  corresponds  to 
weights  optimized  to  minimize  the  energy  in  the  sidelobe  region  of  the  response,  the  left  edge  of 
which  for  this  purpose  is  marked  with  a  O  at  5  dB  above  the  plot’s  bottom  edge.  That  minimization 
is  subject  to  these  several  constraints. 

•  The  boresight  response  must  not  be  less  than  0  dB. 

•  Taper  loss  as  a  line  array  must  not  exceed  0.85  dB. 

•  Response  magnitude  must  not  exceed  —35  dB  anywhere  on  a  grid  of  closely  spaced  frequen¬ 
cies  (or  direction  cosines)  in  the  sidelobe  region  that  extends  rightward  from  the  x  mark 
(plotted  at  that  —  35  dB  level). 

The  first  two  of  these  constraints  are  satisfied  with  equality,  of  course,  and  the  third  is  satisfied  with 
equality  for  some  grid  frequencies. 

This  optimization  specification  results  in  a  response  quite  similar  to  the  Taylor  response  and 
so  indirectly  sheds  light  on  the  Taylor  response  itself:  this  is  more  or  less  the  sense  in  which  the 
Taylor  response  is  roughly  optimal.  Later  we’ll  choose  optimality  criteria  quite  differently,  as  there 
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is  nothing  terribly  special  about  this  combination  of  criteria.  Meanwhile,  here  the  response 

is  for  comparison  purposes  only;  the  optimized  response  is  used  as  the  prototype  for  the  McClellan 
transformation  to  follow. 

McClellan  transformation  ([15])  uses  the  linear-phase  spreading  function  ^3d(  )  X2d(  ) 
of  ([16])  and  is  represented  by  the  contour  plot  on  the  second  line  of  Fig.  [8]  to  map  the  optimized 
ID  response  of  the  top  line  to  a  2D  array  factor  H(X 2d(  ))  with  the  dB  magnitude 

response  plotted  on  two  different  scales  to  the  right  of  the  contour  plot.  In  both  cases  plot  colors 
correspond  to  dB  as  shown  on  the  colorbar  at  the  figure’s  lower  left,  but  the  extent  of  the  center  plot 
in  the  array  plane  matches  the  extent  of  the  contour  plot,  while  the  plot  on  the  right,  though  still 
centered  on  boresight,  covers  an  area  only  about  7.5°  square.  On  the  latter  plot  contours  have  been 
added  at  —1,  —2,  and  —3  dB  to  make  it  easy  to  visually  estimate  typical  beamwidth  measures. 

The  feathery  look  deep  in  the  circular  array-factor  nulls  of  the  rightmost  of  those  plots  is  a 
plotting  artifact,  as  are  the  center  plot’s  apparent  small  discontinuities  at  elevations  near  —35°,  —27°, 
—  12°,  8°,  and  14°.  Similar  discontinuities  are  visible  in  the  rightmost  plot  as  well,  particularly  in  the 
contours.  The  feathers  disappear  when  a  sufficiently  high  plot  resolution  is  used,  but  at  the  expense 
of  course  of  burdensome  plot  computations  and  enormous  plot  files.  The  discontinuities  are  not 
present  in  the  original  plot  but  appear  during  a  file  format  conversion. 

Frequency-warping  function  X^1  (X2d(  ))  shown  using  contours  in  Fig.^is  analogous  to  the 
light  curve  in  Fig.  3j  In  the  latter  the  curve  is  a  function  of  ID  frequency  /,  but  here,  in  Fig.  8j  it  is  a 
function  of  the  2D  array-plane  component  of  3D  spatial-frequency  vector  .  The  level-crossing  dots 
•  on  the  Fig.[3]curve  are  analogous  to  contours  in  Fig.[8]and  convey  exactly  the  same  meaning.  The 
contour  that  marks  FIR  frequency  fS  =  0.4  or  direction  cosine  Xf  =  0.8  lies  along  positions  in  the 
array  plane  where  transformed  array  factor  H  ( X2d  (  ) )  takes  on  the  same  value  that  H(X\v  (fS) ) 
has  at  fS  =  0.4  or  Xf  =  0.8.  Thus  contours  are  level  curves  not  only  for  the  frequency- warping 
function  but  for  the  transformed  array  factor  as  well:  the  path  in  the  array  plane  followed  by  any 
contour  finds  the  transformed  array  factor  holding  absolutely  constant.  This  gives  that  array  factor 
a  very  unusual  look,  with  individual  sidelobes  extending  over  long  distances  in  the  plane. 

The  third  row  of  plots  in  Fig.  [8]  is  like  the  second  except  that  Fourier  pair  ^3d(  )  X2D(  ) 
of  the  second  row  has  been  replaced  with  steered  version  (  )ej2?rsy  X2d((  —  s')  ), 

where  unit  vector  s'  is  at  40°  azimuth  and  30°  elevation.  The  resulting  frequency-warping  function 
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1D1  (X2d  (  ) )  and  array  factor  H (X2d  ( 


) )  respectively  become  steered  functions  X^1  (X2d  ( (  — 
s )  ))  and  i7(X2 d((  —  s')  )).  The  point  made  here  is  simple:  in  those  systems  in  which  the  Mc¬ 
Clellan  transformation  is  used  to  compute  weights  in  real  time,  explicit  phase-shift  steering  of  the 
perhaps  thousands  of  weights  of  taper  h(  )  is  unnecessary.  It  is  enough  to  steer  the  very  few  weights 
of  spreading  function  £30  (  ).  The  McClellan  transformation  handles  the  rest. 

To  put  the  performance  of  the  array  factor  obtained  here  through  the  McClellan  transformation 
into  context,  the  bottom  row  of  plots  in  Fig.  [8]  presents  an  array  factor  for  comparison  for  which 
no  McClellan  transformation  was  used.  Instead  the  corresponding  2D  taper  was  optimized  directly. 
This  array  is  not  hexagonal  in  outline  but  more  or  less  circular:  elements  are  placed  at  all  points  of 
the  X/y/3  spaced  hexagonal  lattice  that  lie  at  distance  A^/l237/3  or  less  from  the  origin.  This  yields 
an  array  with  4507  elements,  making  it  some  3.7%  smaller  than  the  design  of  Fig.  [8]plot  rows  two 
and  three.  The  taper  was  made  invariant  to  rotations  by  multiples  of  60°  and  to  reflections  through 
axes  at  multiples  of  30°  from  the  azimuth  axis.  This  meant  only  404  distinct  weights  required 
optimization,  with  the  twelvefold  symmetry  determining  the  rest. 

The  optimization  maximized  SNR  gain  by  minimizing  noise  gain  while  constraining  the  bore- 
sight  array  factor  to  at  least  unity,  and  additional  constraints  were  imposed  to  bound  the  array-factor 
magnitude  to  no  more  than  —35  dB  outside  an  array-plane  circle  centered  on  the  origin  and  extend- 
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ing  to  2.5°  in  azimuth  and  elevation.  The  latter  required  individual  constraints  at  approximately 
8,000  specific  points  in  the  array  plane,  all  in  (the  closure  of)  a  wedge  covering  1/12  of  an  array- 
factor  period  with  the  rest  of  the  array  plane  taken  care  of  by  the  imposed  symmetries.  The  opti¬ 
mization  specification  was  converted  to  a  second-order  cone  program  (SOCP)  by  our  Opt  matlab 
toolbox  fTO)  and  solved  by  the  SeDuMi  numerical  solver  Q31  12]  (The  latter  tools  were  also  used 
for  the  earlier  ID  optimization.) 

The  optimized  array  factor  in  the  bottom  row  appears  superior  to  the  McClellan-transformation 
design  of  the  second  row  in  several  ways.  It  has  a  very  slightly  narrower  beamwidth.  It  has  simi¬ 
lar  sidelobes  near  the  main  beam  but  considerably  smaller  sidelobes  everywhere  far  from  the  main 
beam.  It  uses  fewer  elements  and  yet  has  some  0.83  dB  more  SNR  gain  by  virtue  of  having  much 
lower  taper  loss.  It  is  not  clear  from  this  one  example  whether  the  McClellan-transformation  ap¬ 
proach  is  seriously  inferior  to  2D  optimization  or  whether  we  just  need  to  be  more  clever  in  applying 
it.  The  next  few  designs  will  move  towards  better  McClellan  designs  in  hopes  of  resolving  this  ques¬ 
tion. 

The  taper-loss  numbers  in  Fig.  [8]  do  demonstrate  nicely  that  taper  loss  is  not  preserved  by  the 
McClellan  transformation  from  ID  to  2D.  From  the  Section [3.2.3  discussion,  the  increase  in  taper 
loss  in  going  from  ID  to  2D  should  be 
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2V2D 1 1  ^2D 1 1 2 
2Vid  ||  ^1D  ||2 


in  terms  of  the  numbers  TVid  and  A^d  of  array  elements  in  ID  and  2D  and  the  sums  ||/md||2  of 
II  ^2D  || 2  of  squared  magnitudes  of  the  ID  and  2D  weights.  Obtaining  zero  taper-loss  increase  would 
require 

1 1  ^2D 1 1 2  =  ^||/,1D||2, 

which  in  this  particular  case  would  amount  to  a  “noise  power  beamwidth”  (in  analogy  to  the  noise 
power  bandwidth  of  a  filter)  of  1.94°.  This  is  fairly  close  to  the  3  dB  beamwidth  obtained  on  the  right 
in  plot  row  two,  so  the  observed  taper-loss  deterioration  in  going  from  ID  to  2D  of  0.96  dB,  a  power 
ratio  of  about  5/4,  likely  arises  largely  from  the  contribution  of  the  large  2D  sidelobe  area  to  ||/i2d||- 
This  suggests  sidelobe  suppression  as  important  to  improve  in  future  McClellan-transformation 
experiments. 


3.3.2  Design  example:  a  larger  spreading  function  allows  beam  spoiling 

To  obtain  more  control  over  Fourier  pair  X2d(  )  and  hence,  through  the  McClellan 

transformation,  over  transformed  array  factor  H(X 2d(  )),  in  Fig. [9] we  increase  the  spatial  extent 
of  ^3d(  )  from  seven  impulses  on  element-location  lattice  points  to  nineteen  impulses  on  two  con¬ 
centric  hexagons  of  points  around  the  origin.  Full  12-fold  symmetry  on  the  impulse  weights  reduces 
the  number  of  design  degrees  of  freedom  to  four,  and  setting  the  overall  scale  and  origin-impulse 
area  to  give  X2B  (  )  a  range  of  [— 1 , 1]  further  reduces  it  to  two.  These  two  are  in  effect  the  ratio  of 

the  side  impulse  area  to  the  corner  impulse  area  on  the  outer  hexagon,  the  “side  ratio”  and  the  ratio 
of  the  corner  area  on  the  outer  hexagon  to  the  corner  area  on  the  inner  hexagon,  the  “corner  ratio’.’ 

Two  frequency- warping  functions  X^  (X2d(  ))  are  presented. 

narrow  beam:  In  the  second  plot  row  the  side  and  corner  ratios  are  each  set  to  1/2  to  give  a  bit 
of  taper  to  the  outermost  impulses.  This  £30  (  )  X2d(  )  pair,  chosen  by  trial  and  error, 
yields  a  narrow  beam  in  transformed  array  factor  H(X 2d(  ))  and  gives  it  a  simple  sidelobe 

structure. 
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wide  beam:  In  the  third  plot  row  the  side  ratio  is  set  to  zero  to  give  pair  £30  (  )  ^  -^2d(  )  a 

“third  band”  structure  m,  which  with  small  sidelobes  will  automatically  force  a  wide,  flat 
beam  in  X2b{  )•  The  corner  ratio  is  set  experimentally  to  about  —0.179  to  obtain  low  outer 
sidelobes  in  X2D  (  ) .  A  wide  beam  in  transformed  array  factor  H ( X2D  (  ) )  results. 

The  top-row  prototype  response  is  again  shown  in  Taylor  and  optimized  versions,  but  now  the 
latter  has  the  taper-loss  bound  loosened  to  1.25  dB  to  give  the  minimization  more  latitude  to  suppress 
sidelobe  energy  in  hopes  of  lowering  the  taper-loss  penalty  in  going  from  ID  to  2D.  This  strategy 
does  not  appear  particularly  successful,  but  further  investigation  would  be  required  in  order  to  say 
anything  definitive  on  this  point.  All  we  can  really  say  here  is  that  the  transformed  2D  sidelobes 
now  appear  at  a  more  generally  reasonable  level  for  the  array  application. 

The  main  purpose  of  this  example  is  to  show  that  a  modest  change  in  the  spreading  function  can 
be  used  to  broaden  or  “spoil”  the  beam  of  transformed  array  factor  H(X 2d(  ))•  Here  the  3  dB 
beamwidth  was  increased  by  a  factor  of  about  3.4.  Unfortunately,  this  spoiling  feature  cannot  be 
used  for  phase-only  transmit  tapers,  because  the  McClellan  transformation  transforms  phase-only 
ID  tapers  to  2D  tapers  that  lack  the  phase-only  property. 

The  taper  loss  shown  for  the  broadened-beam  case  is  over  12  dB,  but  most  of  this  is  simply  an 
unavoidable  geometrical  phenomenon  tied  to  the  beam  broadening  itself  in  a  very  fundamental  way: 
beamwidth  broadening  by  3.4 x  would  naturally  increase  the  portion  of  the  taper  loss  associated 
with  the  main  beam  by  some  20  log10  3.4  «  10.6  dB.  The  taper  loss  here  increases  by  less  than  that 
amount  in  going  from  the  narrow-beam  case  to  the  wide-beam  case,  indicating  that  a  significant  part 
of  the  narrow-beam  taper  loss  comes  from  the  energy  in  the  sidelobe  structure.  Discussion  of  taper 
loss  in  greater  depth  and  an  approach  to  modifying  the  measure  appropriately  for  wide-beam  cases 
is  detailed  in  lfl4l. 

How  much  can  full  optimization  improve  on  the  third-row  design?  The  optimized  design  pre¬ 
sented  in  the  bottom  row  of  Fig.  [^begins  with  array  elements  at  the  4315  points  of  the  X/y/3  spaced 
hexagonal  lattice  that  lie  at  distance  A^/ll91/3  or  less  from  the  origin,  making  the  array  about 
0.4%  smaller  than  the  design  of  Fig.  [9]plot  rows  two  and  three.  The  taper  is  again  constrained  to  be 
twelvefold  symmetric  and  is  optimized  to  minimize  noise  gain  with  additional  constraints  bounding 
the  array-factor  magnitude  to  no  more  than  —35  dB  outside  a  circle  with  a  10° radius. 

Recall  that  the  bottom-row  optimization  in  Fig.  [8]  was  similar  but  additionally  required  the  array 
factor  to  exceed  unity  at  the  origin  so  that  the  minimization  of  noise  gain  would  not  simply  zero 
all  the  weights.  Here  the  array  factor  is  instead  required,  everywhere  inside  a  circle  of  radius  3.3°, 
approximating  the  3  dB  beamwidth,  to  exceed  the  broadened  third-row  array  factor  but  to  not  exceed 
unity.  Sandwiching  the  main  beam  between  these  bounds  forces  a  beam  shape  very  similar  to  that 
shown  in  the  third  row.  We  can  then  compare  the  sidelobe  levels  to  get  a  sense  of  how  much  better 
true  optimization  is  than  the  McClellan  approach  used  for  row  three. 

The  optimized  fourth-row  taper  turns  out  to  be  quite  comparable  to  the  third-row  McClellan- 
transformed  design  in  3  dB  beamwidth,  in  taper  loss,  and  in  first-null  beamwidth,  but  its  sidelobe 
performance  is  so  much  better  than  the  McClellan-transformed  design  that  a  different  color  scale 
was  used  to  plot  it,  one  that  covers  50%  more  dB. 

So  while  the  ready  availability  of  simple  beamspoiling  is  an  attractive  feature  of  the  McClellan 
transformation,  the  question  of  ultimately  obtainable  performance  for  a  given  number  of  elements 
appears  to  be  serious  and  suggests  the  next  experiment. 

3.3.3  Design  example:  trading  ripple  for  beamwidth 
The  Fig.  [TO] design  experiment  has  four  aims. 
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1 .  Make  the  array  shape  more  circular  by  mov¬ 
ing  away  from  hexagonally  shaped  spread¬ 
ing  functions. 

2.  Experiment  with  a  less  enormous  array. 

3.  Examine  the  tradeoff  between  ripple  and 
beamwidth  in  the  design  of  spreading  func¬ 
tions. 

4.  Improve  the  design  process  by  jointly  opti¬ 
mizing  the  ID  prototype  filter  and  the  2D 
McClellan-transformed  array  factor. 

Aim  [I]  is  addressed  through  constructing  a 
spreading  function  )  =  x^{  )  with  impulses 
located  on  the  31  lattice  points  shown  in  Fig.  [7] 
These  points  form  three  concentric  hexagons,  all 
centered  on  the  origin,  with  sides  of  two,  three, 
and  four  points,  but  with  the  six  corner  points 
of  the  outer  hexagon  removed.  This  is  the  least 
number  of  impulses  that  yields  nonhexagonal  ar¬ 
ray  shapes  while  retaining  twelvefold  symmetry. 
The  basis  functions  t 1  (  ) ,  £2  (  ) , . . .  and  therefore 
of  the  array  itself  have  dodecagonal  (12- sided)  out¬ 
lines.  They  are  hexagons  in  the  spirit  of  Fig.  [5]but 
with  the  six  corners  “filed  down’.’  The  hope  was 
that  by  not  wasting  elements  on  the  corners  where 
not  really  needed  the  2D  taper  losses  could  be  im¬ 
proved. 

Aim  [2]  is  addressed  by  setting  M  —  10,  result¬ 
ing  in  an  array  size  equal  to  the  number  of  im¬ 
pulses  in  f  10  (  )  or  2461  elements. 

Aim  [3]  is  addressed  in  the  Fig.  [TO]  experiment 
by  carrying  out  two  versions  of  the  design.  The 
spreading-function  pairs  x^(  )  ^  X2d(  )  of 
plot  rows  two  and  three  are  each  optimized  to  have 
X2D  (0)  =  1  and  for  further  from  the  origin 
than  some  stopband  radius  to  have  X2d(  )  con¬ 
tained  in  some  interval  [—1,6]  with  bound  b  mini¬ 
mized.  It’s  a  shifted  and  scaled  version  of  a  classic 
lowpass-filter  design  problem,  here  with  a  single¬ 
point  passband  at  the  origin.  The  difference  is  that 
here  the  usual  near-zero  stopband  is  replaced  by  a 
region  near  —1  for  which  the  stopband  terminol¬ 
ogy  we  will  use  is  at  best  appropriate  by  analogy. 

There  is  a  tradeoff  between  the  stopband- 
radius  design  parameter  and  the  optimal  value  of  b 
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that  results.  The  small  radius  of  Fig.  10  row  two  yields  a  relatively  high  bound  b  «  cos(27tx0.1764) 
so  that  the  stopband-ripple  swings  are  larger.  The  high  radius  of  row  three  yields  a  very  low  bound 
b  «  cos(27t  x  0.4666),  making  the  stopband-ripple  swings  very  small.  In  either  case  the  obtained 
bounds  b  on  X2d(  )  imply  bounds  b\  here  0.1764  and  0.4666  in  the  two  cases,  on  frequency¬ 
warping  function  X1d1(X2d(  )),  and  the  latter  bounds  are  then  used  to  set  up  the  top-row  opti¬ 

mization  of  the  ID  prototype  filter. 

Aim[4|  joint  optimization  of  the  prototype  filter  and  the  McClellan-transformed  array  factor, 
is  accomplished  by  optimizing  the  common  set  of  weights  that  determine  them  both  using  some 
specifications  in  ID  and  some  specifications  in  2D. 


•  The  prototype  filter’s  DC  gain  H(X id(0)),  which  is  also  the  boresight  array  factor  H(X 2d(0)), 
must  be  at  least  unity  or  0  dB, 

•  prototype  filter  H(Xm(fS))  is  limited  to  —30  dB  stopband  ripple  for  0.08  <  fS  <  b\ 

•  prototype  filter  H(X\v(fS))  is  limited  to  —50  dB  stopband  ripple  for  b'  <  fS  <  0.5,  and 

•  the  2D  taper  loss  of  array  factor  H(X 2d(  ))  is  minimized. 


The  idea  is  that  the  McClellan  transformation  guarantees  that  the  ripple  magnitudes  are  the  same  in 
ID  and  2D,  but  bounding  that  magnitude  is  easier  in  ID  because  far  fewer  constraints  are  needed 
there.  The  taper  loss  that  matters,  however,  is  the  2D  taper  loss. 

Joint  ID  and  2D  optimization  would  seem  to  be  ideal,  but  embedded  in  that  design  is  the  con¬ 
straining  effect  of  the  hand-picked  spreading-function  stopband  radius.  By  determining  the  spread¬ 
ing  function,  that  radius  determines  the  basis  set  used  to  build  the  2D  array  factor  and  therefore 
the  subspace  to  which  its  design  is  restricted.  Aim [3] is  therefore  key  here:  is  it  better  for  the  ID 
prototype  to  have  a  smaller  outer  stopband  to  get  the  narrowest  main  lobe  possible?  Or  would  a 
larger  outer  stopband  and  a  slightly  wider  main  lobe  be  preferable? 

Comparing  Fig.  [10]  rows  two  and  three  answers  this  quite  clearly  for  the  design  of  conven¬ 
tional  pencil  beams.  A  small  spreading-function  stopband  radius  is  preferable,  even  though  a  large 
spreading-function  ripple  results.  That  large  ripple  requires  a  concomitantly  low  outer- stopband 
cutoff  in  the  ID  prototype  filter,  which  leads  to  a  larger  ID  main  beam,  but  the  latter  penalty  is  more 
than  overcome  when  the  McClellan  transformation  maps  that  main  beam  into  a  much  smaller  2D 
beam.  Beamwidth  is  reduced  in  the  favored  approach  by  over  20%,  and  the  width  of  the  —30  dB 
inner  sidelobe  shelf  is  markedly  reduced  as  well. 

The  optimized  design  of  the  bottom  row  of  Fig.  10  is  also  twelvefold  symmetric  with  —30  dB 
inner  sidelobes,  with  —50  dB  outer  sidelobes,  and  with  the  transition  between  constraint  regions 
occurring  at  4°  and  8°  radii  for  comparability  to  the  design  of  row  two.  Relative  to  that  McClellan- 
transformed  design,  however,  the  optimized  bottom-row  design  has  a  very  slightly  narrower  beam, 
clearly  lower  average  sidelobe  levels,  and  0.37  dB  more  SNR  gain,  and  it  achieves  this  with  over 
14%  fewer  elements,  as  it  has  elements  only  at  the  2107  points  of  the  A/a/3  spaced  hexagonal  lattice 
that  lie  at  distance  Ay/579/3  or  less  from  the  origin.  The  taper  is  optimized  to  minimize  noise  gain 
subject  to  constraints  on  peak  sidelobe  levels. 


3.3.4  Design  example:  small  spreading  function  vs.  large 

If  questions  of  array  shape  or  beam  spoiling,  etc.  are  not  at  issue,  is  there  any  inherent  performance 
advantage,  in  the  beamshaping  application,  in  a  spreading  function  comprising  a  smaller  or  larger 
number  of  impulses?  The  design  examples  of  Figure  |TT]  address  this  question.  The  hexagonal 
supports  of  the  spreading  functions  used  in  the  designs  of  the  second  and  third  rows  contain  seven 
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and  37  points  respectively.  These  are  used  with  M  —  30  and  M  =  10  respectively  so  that  by  starting 
with  prototype  filters  of  different  lengths,  the  McClellan  transformation  yields  a  hexagonally  shaped 
array  of  2791  elements  in  each  case  in  spite  of  the  difference  in  spreading-function  size.  Each 
spreading  function  is  optimally  designed  exactly  as  in  the  Fig.  [Redesigns  discussed  in  Section [J. 3. 4 
except  that  in  each  Figure  |TT]  design  the  spreading-function  stopband  radius  was  set  to  10°.  The 
very  different  spreading-function  sizes  then  naturally  led  to  very  different  ripple  ranges  [— 1 ,6]  for 
spreading  function  X2d(  )  in  the  two  cases,  and  therefore  to  very  different  stopbands  [b\  0.5]  for 
frequency-warping  function  X1d1(X2d(  ))  as  well.  The  b' values  of  0.0817  and  0.2361  shown  in 

Figure  [IT]  for  the  second-  and  third-row  designs  then  become  the  stopband  cutoff  frequencies  used 
to  optimize  the  associated  prototype  filters  of  the  first  row  of  Figure  [11] 

A  closer  look  at  the  significance  of  the  10°  spreading-function  stopband  radius  is  enlightening. 
That  radius  actually  divides  the  plane  into  an  inner  region  where  the  [0,  b ']  interval  in  ID  frequency 
for  the  prototype  filter  maps  monotonically  to  radius  in  2D  and  an  outer  region  where  there  is  no 
such  monotonicity,  where  ID  frequency  moves  back  and  forth  across  the  [&'  0.5]  stopband  in  various 
ways  as  2D  radius  increases  beyond  10°.  Because  of  that  complexity  in  the  outer  region,  it  makes 
little  sense  to  design  for  any  particular  detailed  stopband  behavior  there,  and  these  designs  simply 
use  an  upper  bound  on  array-factor  magnitude  there.  The  monotonicity  in  the  inner  region,  however, 
makes  it  reasonable  to  attempt  some  shaping  of  the  2D  stopband  behavior  there  through  shaping  of 
the  ID  prototype  filter’s  stopband  behavior  in  the  associated  interval. 

In  the  Fig.  [TO] designs  discussed  earlier  in  Section  3.3.4  the  inner  sidelobes  were  set  to  a  con¬ 
stant  level  of  —30  dB,  while  the  outer  sidelobes  were  bounded  to  —50  dB  and  below.  Here  in 
Figure  [IT]  however,  only  the  innermost  sidelobe  is  permitted  to  reach  —  30  dB,  and  in  the  first-row 
ID  prototypes  the  bounding  level  used  to  control  the  other  inner  sidelobes  is  ramped  linearly  down 
to  the  —50  dB  outer- sidelobe  level.  The  resulting  decreasing  heights  of  the  inner  sidelobe  rings  in 
the  McClellan-transformed  array  factors  is  a  generally  desirable  feature,  and  this  experiment  also 
satisfies  some  basic  curiosity  as  to  whether  such  a  feature  might  reduce  the  penalty  suffered  by 
the  McClellan  approach  in  taper  loss  and  SNR  gain  relative  to  optimized  designs.  (The  answer: 
perhaps,  but  not  dramatically.) 

The  bottom-row  array  factor  was  for  the  same  hexagonal  array  and  was  directly  optimized 

•  to  lower  bound  its  boresight  magnitude  at  0  dB, 


to  upper  bound  its  magnitude  at  —30  dB  from  3.5°  to  4.1°  in  radius, 


•  to  upper  bound  its  magnitude  with  a  decibel-scale  ramp  extending  from  —30  dB  at  4.1°radius 
to  —50  dB  at  10°  radius, 


•  to  upper  bound  its  magnitude  at  —50  dB  everywhere  outside  10°  radius,  and 

•  to  minimize  taper  loss. 

The  differences  in  taper  loss  among  the  three  designs  appears  entirely  accounted  for  by  their 
differences  in  beamwidth.  The  taper  loss  is,  as  expected,  best  for  the  fully  optimized  design,  by  a 
substantial  margin  of  0.72  dB.  That  fully  optimized  design  also  had  far  lower  sidelobe  levels  over 
well  over  half  the  hemisphere,  though  such  levels  are  quite  pointless  because  hardware  weighting 
errors  would  prevent  their  realization  in  practice.  Better  would  be  to  specify  the  optimization  in  a 
way  that  used  the  extra  degrees  of  design  freedom  to  more  useful  effect. 

The  two  McClellan-transformed  designs  were  surprisingly  close  in  performance,  differing  in 
taper  loss  by  only  0.14  dB.  It  was  no  surprise,  however,  that  the  slightly  better  of  the  two  was 
the  design  based  on  the  small  spreading  function,  simply  because  that  optimization  worked  with 
M+l  =  31  basis  functions  as  opposed  to  11. 
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3.3.5  Design  example:  a  tilted-ellipse  beamshape 

Reddy  and  Hazra  lfl5l  used  a  McClellan  transformation  to  create  elliptical  passbands  or  beamshapes 
with  the  ellipses  oriented  arbitrarily  in  the  plane.  They  assumed  a  square  sampling  (element  layout) 
lattice  in  the  plane,  however,  and  they  considered  only  a  very  simple  nine-sample  spreading  func¬ 
tion.  Let  us  see  what  we  can  do  along  these  lines  but  for  the  hexagonal  lattice  and  using  optimization 
tools  to  favor  high-performance  beamshapes. 

The  design  experiment  presented  in  Fig.  [I2|  investigates  broadening  the  beam  in  one  arbitrarily 
chosen  direction  only.  To  do  this  it  combines  the  geometry  of  the  Fig.  [10] design,  in  particular 

•  the  31 -point  spreading-function  support  and  the  resulting  dodecagonal  array  shape  and 


•  an  array  sized  by  setting  M  — 10  to  yield  2461  array  elements, 
with  the  basic  structure  of  the  Fig.  [IT] design,  in  particular 

•  inner- sidelobe  heights  that  ramp  down  from  —30  dB  at  a  fixed  ID  cutoff  frequency  of  fS  = 
0.098  to  —50  dB  at  the  outer- sidelobe  bandedge  b'  that  results  from  the  spreading-function 


optimization  as  discussed  in  Section  3.3.3 


These  features  from  prior  designs  are  combined  with  one  new  aspect: 


•  the  stopband-edge  circle  of  10°  radius  used  to  constrain  previous  spreading-function  designs 
is  here  replaced  with  a  stopband  edge  that  is  elliptical  in  the  two  dimensions  of  .  That 
ellipse 


-  has  a  minor  radius  equal  to  the  10°  radius  of  the  circle  used  in  earlier  designs, 

-  has  a  major  radius  double  that — the  doubling  is  in  ||  ||  rather  than  angle  so  that  the 

major  radius  is  slightly  larger  than  20°  in  angle — and 

-  the  ellipse  is  placed  in  the  plane  with  its  semi-major  axis  at  an  angle  to  the  hori¬ 
zontal  (azimuth)  axis  that  begins  at  zero  in  row  two  of  Fig.  [12]  and  that  increments  in  5° 
increments  in  successive  rows  of  the  figure. 

Obtaining  angles  at  5°  increments  beyond  30°  is  a  matter  of  taking  one  of  the  designs  in  the  [0°,  30°] 
interval  and  either  rotating  its  spreading-function  weights  by  a  multiple  of  60°  or  reflecting  those 
weights  through  an  axis  at  some  multiple  of  30°  from  horizontal.  Establishing  the  behavior  of  the 
design  over  seven  angles  at  5° increments  in  the  [0°,  30°]  interval  therefore  effectively  establishes  its 
behavior  at  72  multiples  of  5°  spaced  evenly  around  the  circle. 

The  first  row  of  Figure [l2]presents  the  overlaid  magnitude  responses  of  the  separately  optimized 
ID  prototype  filters  of  all  seven  designs,  along  with  their  stopband  bounding  lines.  The  seven 
prototypes  are  sufficiently  similar  that  it  appears  entirely  practical  to  choose  one  fixed  prototype 
design  to  use  at  all  angles. 

Figure[l2]breaks  from  this  report’s  pattern  of  showing  the  entire  hemisphere  in  the  center  column 
and  instead  covers  only  about  ±20°  relative  to  boresight  in  each  direction,  primarily  to  simplify 
examination  of  the  stopband-boundary  ellipse,  shown  in  red,  and  the  inner- sidelobe  region.  As  in 
earlier  designs  the  outer  sidelobes  in  the  rest  of  the  hexagonal  period  of  array  factor  X2D  (  )  simply 
comprise  extended  ripple  ridges  that  peak  at  —50  dB  as  before.  Those  outer  sidelobes  move  and 
change  their  ridge  geometries  with  ellipse  angle  but  are  otherwise  not  terribly  interesting. 

The  third  column  of  Figure  [12]  shows  the  center  portion  of  the  main  beam  up  very  close,  closer 
than  for  other  designs,  so  that  the  shapes  of  the  contours  can  be  visually  compared  to  the  shapes  of 
the  stopband-boundary  ellipses  in  the  center  column. 
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As  the  angle  of  the  stopband-boundary  ellipse  is  incremented,  what  can  be  seen  is  that  even 
though  the  —  3  dB  main-beam  contour  is  not  explicitly  shape  controlled  in  the  optimization  and 
obtains  its  shape  only  indirectly  from  the  stopband-boundary  ellipse,  that  main-beam  contour  gen¬ 
erally  becomes  elliptical  and  tilts  right  along  with  the  stopband-boundary  ellipse.  However,  that 
angle  tracking  is  neither  perfect  nor  does  it  maintain  the  2  : 1  axis  ratio  consistently.  The  20°  and 
25°  tilts  in  particular  have  main-beam  ellipses  that  are  visibly  less  eccentric.  The  designs  for  each 
of  these  two  tilts,  and  to  a  lesser  extent  the  design  for  the  15°  tilt  as  well,  also  show  changes  in 
how  the  center-column  stopband-boundary  ellipse  relates  to  the  sidelobe  ridges  immediately  around 
it.  Somehow  the  stopband-boundary  ellipse  “loses  traction”  on  the  design  a  bit  at  these  middle  tilt 
angles. 

Roughly  3  dB  of  the  taper  losses  shown  in  the  third  column  of  course  are  the  organic  result 
of  realizing  the  intended  beam  broadening  by  a  factor  of  two  in  area  and  do  not  represent  some 
sort  of  greatly  increased  design  inferiority.  The  last  two  tilts,  25°  and  30°,  result  in  very  different 
mainbeam  ellipticities  and  correspondingly  very  different  taper  losses.  This  shape  sensitivity  makes 
these  numbers  much  less  useful  here  than  in  earlier,  pencil-beam  designs.  They  are  included  here 
for  completeness. 

The  take-away  lesson  from  the  designs  for  this  series  of  seven  tilts  is  that  the  general  idea  of 
simply  changing  the  spreading-function  weights  to  effect  beam  broadening,  whether  in  ID  as  here 
or  in  2D  as  in  the  Fig.[9]design,  appears  to  be  quite  viable.  Further,  it  appears  that  given  a  particular 
desired  elliptical  beam  shape,  a  single  prototype  design  is  likely  viable  for  all  ellipse  tilt  angles.  The 
detailed  design  procedure  used  here  is  probably  not  quite  adequate  in  the  elliptical-beam  case,  but 
it  should  not  be  particularly  difficult  to  modify  the  optimization  formulation  to  better  preserve  the 
beam  shape  as  the  ellipse  tilt  is  varied,  but  we  can  safely  leave  demonstration  of  that  to  the  relevant 
application  design. 

In  a  real  application  making  the  tilt  a  continuously  variable  parameter  is  apt  to  be  an  attractive 
option,  and  the  methods  recently  proposed  by  Yeung  and  Chan  [16 1  and  by  Shyu  et  al.  [17]  to  com¬ 
bine  the  McClellan-transformation  approach  to  2D  filters  with  the  Farrow  structure,  or  variations 
thereon,  for  realizing  variable  FIR  filters  should  be  examined  closely. 


4  Conclusions 


Broadly  speaking,  the  McClellan  transformation  seems  well  suited  to  quick  design  of  2D  array  ta¬ 
pers  of  medium  to  large  size.  In  this  report  most  of  the  ID  prototype  filters  transformed  are  designed 
using  optimization  approaches,  but  the  real  purpose  of  such  optimization  here  is  to  learn  something 
about  obtainable  limits  to  performance.  In  practice  the  strength  of  McClellan  transformation  is  that 
it  permits  easy  creation  of  large  2D  tapers  by  using  ordinary  design  tools  to  design  a  ID  taper  and 
using  ordinary  2D-filter  design  tools  to  design  a  very  small  2D  spreading  function.  Reasonable  re¬ 
sults  can  even  be  obtained  using  very  small  spreading  functions  designed  by  hand  and  simple  Taylor 
or  Remez-exchange  design  of  prototypes. 

This  report’s  design  experiments  reveal  only  one  fundamentally  important  limitation  of  the 
McClellan-transformation  approach  to  2D  taper  design:  taper  loss.  Every  combination  of  spreading- 
function  and  prototype-filter  design  criteria  yielded  a  taper  loss  at  least  0.7  dB  greater  than  obtained 
for  the  same  size  array  through  direct  optimization  of  the  taperj^]  This  is  a  significant  limitation 
and  may  limit  the  use  of  the  McClellan  transformation  in  final  designs  to  those  rare  settings  in 
which  the  approach’s  enabling  of  real-time  changes  to  the  2D  taper  through  simple  modifications  to 


2 In  the  Fig.  llOl design  of  Section 


3.3.3 


SNR  performance  only  0.37  dB  worse  than  the  fully  optimized  reference 


design  was  achieved,  but  the  latter  design  used  some  14.4%  fewer  elements,  which  amounts  to  another  0.67  dB. 
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the  spreading  function  is  of  some  special  benefit.  Also,  in  preliminary  design  work  in  which  opti¬ 
mization  code  capable  of  handling  large  2D  problems  is  unavailable,  the  McClellan-transformation 
approach  might  reasonably  be  used  to  get  early  estimates  of  pattern  shape,  perhaps  including  in  SNR 
or  beamwidth  estimates  some  rough  fudge  factor  to  compensate  for  further  improvements  likely  to 
be  obtained  through  eventual  full-out  optimization. 

Mersereau  and  Mecklenbrauker  J5|  assert  that  the  McClellan  approach  can  save  an  enormous 
number  of  coefficient  multiplies  in  implementation.  This  claim  needs  to  be  interpreted  carefully  in 
the  array  context,  however.  Basically,  the  defining  recursion  for  the  Chebyshev  polynomials  leads 
to  a  structurally  similar  recursion  at  the  heart  of  the  McClellan  transformation.  The  latter  recursior0 
can  take  any  of  three  related  forms,  each  with  an  apparently  natural  application  to  arrays. 

Receive-array  implementation:  An  FIR-filter  recursion  in  Chebyshev  direct  form,  as  presented  in 
Fig.  [I]  can  incorporate  the  McClellan  recursion  into  the  real-time  computation  of  a  receive- 
array  beam  sum.  Each  prototype-filter  coefficient  is  used  in  multiplication  just  once. 

Transmit-array  implementation:  An  FIR-filter  recursion  in  reversed  Chebyshev  form,  as  pre¬ 
sented  in  Fig.  [2j  can  incorporate  the  McClellan  recursion  into  the  real-time  computation  of 
element-drive  signals.  Each  prototype-filter  coefficient  is  used  in  multiplication  just  once. 

Offline  design,  for  transmit  or  receive:  A  McClellan  recursion  can  be  used  to  compute  the  taper 
weights  offline,  and  those  weights  can  then  be  used  conventionally  to  compute  receive  beam 
sum  or  transmit  element  drives  in  real  time. 

The  relatively  obvious  third  approach  treats  the  McClellan  transformation  as  an  easily  computed 
design  procedure  unrelated  to  issues  of  implementation.  The  first  two  approaches  build  the  Mc¬ 
Clellan  recursion  right  into  the  implementation  so  that  it  operates  on  actual  signals  and  not  just  on 
taper  weights,  so  the  taper  weights  are  never  explicitly  computed  at  all.  This  saves  many  coeffi¬ 
cient  multiplies  if  one  counts  only  multiplications  by  the  weights  of  the  prototype  filter.  In  the  first 
two  approaches  above,  however,  the  2D  spatial  convolutions  with  the  spreading  function,  denoted 
by  x  in  Figs.  |T| and [2] but  by  x^(  )  later  in  array- specific  discussion,  dominate  the  required  im¬ 
plementation.  Some  careful  application-specific  operation  counting  would  seem  necessary  before 
pronouncing  on  any  implementation  savings  to  be  had. 
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Fig.  8  Spreading  0dB 
function  (2nd  row 
left)  maps  ID  re- 
sponse  or  taper 
(top  row)  into  2D  dB 
(2nd  row  center  & 
right).  Steered  -30  dB 
spreading  steers 
array  (3rd  row).  -40  dB 
Compare  to  taper 
directly  optimized  -50  dB 
in  2D  (4th  row). 

See  SectionlaBTI  _60  dB 
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Fig.  10  A  spread¬ 
ing  function  of  31  el¬ 
ements  and  a  proto¬ 
type  filter  of  length 
31  are  optimized  two 
ways:  with  tiny  and 
huge  stopband  radii. 
An  optimal  planar  ar¬ 
ray  is  included  for 
comparison.  See 
Section  13331 
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Fig.  11  Array  of 
2791  elts  as  hex 
of  max  width  61 
elts.  Spreading 
fns  of  seven  and 
37  elts  (2nd,  3rd 
rows)  and  an  op¬ 
timal  design  for 
reference  (right). 
See  Sec.  13.3.41 


-15  dB  -10  dB  - 
-30  dB  -20  dB  - 
-45  dB  -30  dB  - 
-60  dB  -40  dB  - 
-75  dB  -50  dB  - 


30 


filter  of  length  21,  <  —  50  dB  outside  of  2:1  ellipse 


Fig.  12  (2  pages)  The  Fig.[To]row-three  McClellan 
geometry  maps  a  length-21  prototype  to  a  2461- 
element  dodecagonal  array,  here  usingan  ellipti¬ 
cal  stopband  boundary.  See  Section [3~3~5l 
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—  1,  —2,  and  —3  dB  contours 
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